CARS Spectral Fitting with Multiple Resonant Species 

Using Sparse Libraries 


Andrew. D. Cutler* and Gaetano Magnotti 1 * 
The George Washington University, Mechanical and Aerospace Engineering, 

1 Old Oyster Point Road, Suite 200, Newport News, VA, 23602 

Abstract 

The dual pump CARS technique is often used in the study of turbulent flames. Fast and accurate algorithms are 
needed for fitting dual-pump CARS spectra for temperature and multiple chemical species. This paper describes the 
development of such an algorithm. The algorithm employs sparse libraries, whose size grows much more slowly 
with number of species than a conventional library. The method was demonstrated by fitting synthetic 
“experimental” spectra containing 4 resonant species (N 2 , 0 2 , H 2 and C0 2 ), both with noise and without it, and by 
fitting experimental spectra from a H 2 -air flame produced by a Hencken burner. In both studies, weighted least 
squares fitting of signal, as opposed to least squares fitting signal or square-root signal, was shown to produce the 
least random error and minimize bias error in the fitted parameters. 

I. Introduction 

Coherent anti-Stokes Raman spectroscopy (CARS) is a non-linear spectroscopic technique in which three laser 
beams, two pump beams and a Stokes beam, are focused and crossed at their focal point, and a fourth signal laser 
beam is generated at the intersection through a four-wave mixing process. Resonances associated with Raman 
active molecular rotational-vibrational transitions strengthen the CARS signal in a manner dependent on gas 
composition and temperature. The wavelength of excitation or “Raman shift” corresponds to the difference in 
energy between a pump beam and a Stokes beam (longer wavelength) photon. In the broadband CARS technique 
one of the laser beams, the Stokes beam, is spectrally broad whereas the other pump beams are spectrally narrow, 
and this enables multiple resonances to be excited. The signal beam is spectrally broad and carries with it the 
spectral signature of the gases present (over some range of wavelengths). This signal beam is dispersed by a 
spectrometer and spectra are compared to analytical models to obtain information on composition and temperature. 
In the dual-pump CARS technique, originally developed by Robert Lucht and coworkers, 1 2 the two pump beams are 
at different wavelengths and two different ranges of Raman shift may be accessed. 

In recent years, efforts have been made to utilize the CARS technique to obtain experimental data for 
development and validation of computational models of supersonic combustion. Data sets have been acquired in a 
hydrogen-fueled supersonic combustor using the CARS technique 3 to measure temperature, and the dual-pump 
CARS technique has been used to measure temperature and composition; 4 in Ref. 3 the resonant species was N 2 
and in Ref. 4 they were N 2 , 0 2 , and H 2 . Both mean flow and turbulence statistics (variances and covariances) were 
derived from the data, although the uncertainty in the latter was high due both to instrument error and to the small 
number of measurements on which to base the statistics. More recently, data have been acquired in a hydrogen- 
combustion heated supersonic coaxial jet with the intent to obtain turbulence statistics of lower uncertainty. 5 These 
efforts required analysis of large numbers of experimental spectra containing up to three different species. 

Analysis of dual-pump CARS experimental spectra is typically performed with the Sandia CARSFT code as 
modified by Lucht 1 2 6 and O’ Byrne. 4 The CARSFT code is able to accurately compute spectra for a range of 
molecular species and, in the fitting mode, to read an experimental spectrum and fit it to a theoretical spectrum by 
repeatedly computing spectra and minimizing the residual between the theoretical and experimental spectra. This 
fitting process becomes very slow when the number of molecular species increases, due to a large spectral range, 
many energy level transitions, and a large number of iterations. Additionally, when fitting multiple variables, 
CARSFT has problems in finding a global minimum, sometimes returning a result dependent upon the initial guess. 
At NASA Langley, with the exception of Tedder et al., 5 CARS data has usually been fitted using libraries of 
spectra/ 4 Libraries are computed ahead of time using CARSFT. Theoretical spectra are then interpolated from the 
library rather than being calculated from scratch. However, conventional library size scales to the power of the 
number of parameters varied (temperature plus resonant chemical species) and tends to become uneconomically 
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large for more than three fitted species. Recently, dual-pump broadband CARS techniques have been developed 
that employ broadband dye lasers of larger spectral width and simultaneous detection of up to 6 resonant species 
was demonstrated; 7 such spectra could not be fit using conventional library methods and are difficult if not 
impossible to fit using CARSFT in its fitting mode. 

Recently, we have developed fitting algorithms that use sparse libraries (as opposed to the “fully-packed” 
libraries used previously) to enable fitting of a larger number of resonant chemical species and rapid fitting of large 
experimental data bases acquired in combustion turbulence studies. The objective of the current work is to describe 
these techniques and their application to fitting dual-pump CARS spectra containing N 2 , 0 2 , H 2 and C0 2 . 

II. Method 

The fitting method employs a library of theoretical spectra calculated separately, using CARSFT or some other 
similar code. Theoretical spectra are interpolated from the library at each iteration of the fitting algorithm. The 
algorithm fmds the parameters of the fit, which are the temperature, mole fractions of resonant species, and a 
horizontal wavelength shift that minimize a residual function between the experimental and theoretical spectra. 

A Library Structure 

Previously, libraries used in fitting CARS spectra were fully populated on a regular grid of temperature and 
species mole fractions, with size (number of spectra) that scaled with the power of the number of resonant chemical 
species (e.g.. Ref. 4). The new method takes advantage of the fact that the resonant spectra of the different species 
generally do not overlap very much. (Indeed, when doing dual-pump CARS the second pump frequency offset is 
typically chosen so as to minimize overlap.) The number of points in the sparse grid is effectively reduced 
compared to the fully populated grid, while interpolation from the sparse library is of comparable accuracy to 
interpolation from the fully packed library, except at wavelengths where resonances of multiple species overlap. 
The library is only sparse in the species concentration coordinate directions; that is to say, the sparse grid of species 
concentrations is fully repeated at each library temperature. 



Figure 1. (a) Library grid points in weighted coordinates in the N 2 and 0 2 plane (10 points and 3 nodes along 
each line); (b) typical spectra along an N 2 line of the library at compositions indicated by the red points in (a), 
normalized by non-resonant background susceptibility 

The species concentration coordinates are chosen so that the variation of a given species affects primarily the 
portion of the CARS spectrum that is resonant for that species, leaving the rest of the spectrum relatively unchanged. 
Thus, the grid is orthogonal in CARS non-resonant susceptibility weighted fraction coordinates defined below: 

y n — ^n%NR,n / % NR 

Here, x n is the mole fraction of species w, N d is the number of resonant species, species N d + 1 is the non-resonant 

species (thus, V * Nd+l x — 1 X Xnryi is the CARS non-resonant susceptibility of species /?, and %nr is the non-resonant 
»= i n 

susceptibility of the mixture as defined below: 
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_ ^ Nd + 1 

XNR ~ 2^1 n=l X n%NR,n 

The location of grid points in a plane of a typical sparse grid in which N 2 and 0 2 vary, plotted in weighted 
fraction coordinates, is shown in Figure 1(a). For each species there is a line of points at the “nodes” of the other 
species; for the case shown in the figure there are three nodes for each species and 10 points along each line giving a 
total of 60 points. By comparison, in the comparable hilly packed library each point is a node and there would be 
10x10 = 100 points for the example of Figure 1(a). It is assumed that there always exists an additional species that 
has no resonances in the range of wavelengths considered, and that makes up the sum of the mole fractions of all the 
species to 1; in this example the non-resonant species is water vapor. Some points in the library define compositions 
that have negative amounts of the non-resonant species. Spectra at these physically unrealizable compositions can 
none-the-less be computed and may be required for interpolation to realizable compositions in instances when the 
composition contains only a small amount of the non-resonant species. 

The library spectra signal intensities are normalized by the non-resonant susceptibility of the mixture, Xnr , prior 
to interpolations. Figure 1(b) contains plots of the square root of the normalized signal intensity with N 2 , 0 2 , H 2 , and 
C0 2 resonant (see Section El for more details of this calculation) along a N 2 coordinate line in the library (the red 
circles) shown in Figure 1(a). In these spectra the shape varies only in the region where N 2 is resonant, which would 
not be the case if, say, the mole fraction of N 2 was varied at constant mole fraction of the other resonant species. 

Table 1 Library size at each temperature given 10 points and 3 nodes along each line. 


no. 

chem 

species 

sparse 

library 

size 

packed 

library 

size 

ratio 

packed to 
sparse 

1 

10 

10 

1.00 

2 

60 

100 

1.67 

3 

270 

1000 

3.70 

4 

1080 

10000 

9.26 

5 

4050 

100000 

24.69 


Library sizes (at each temperature level) for typical fully-packed libraries with 10 points per coordinate direction 
are compared to the sizes of a sparse library with 10 points per line and 3 nodes in Table 1. Computational times of 
fitting algorithms involving libraries scale roughly with library size, as does computer memory usage. The 
computational savings of using a sparse library become substantial for three chemical species or more, and make the 
fitting of up to at least 5 chemical species practical (which is impractical with a fully-packed library). 


B 


Interpolation from the Library 
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Figure 2. Weight factor distributions derived from the library and used in interpolation of spectra from the 

library 

Before spectra can be interpolated from the sparse library, a weight function distribution is determined for each 
resonant species, W n . This function is an estimate of the typical contribution to the intensity of the spectra of a gas 
mixture from that particular species. For each resonant species the average normalized spectrum is formed by 
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averaging (across all temperatures) the library spectra that have the maximum mole fraction of that species and zero 
of the other resonant species. The weight functions are the average spectra renormalized at each wavelength by the 
sum of the average spectra for all resonant species; thus the sum of the weight functions of all the species is 1. 
Typical weight function distributions computed for the library with N 2 , 0 2 , H 2 , and C0 2 resonant are shown in 
Figure 2. The origins of the various peaks and valleys may be identified by comparing to the typical spectra in Fig. 
1(b), where the resonances in the various species can be identified. Thus, for example, the weight function for N 2 
tends to 1 at wavelengths where N 2 is resonant, has some sharp dips where H 2 resonances overlap the N 2 resonances, 
tends to zero where other species are resonant and N 2 is not resonant, and where no species are resonant has the 
same value (-0.25) for all species. 

The method for interpolating a spectrum at a particular library temperature to a desired composition has several 
steps. A complete spectrum (including all the resonances), T m is interpolated to the desired composition from all the 
library lines in a given species coordinate direction, one interpolated spectrum for each species. Thus, in Fig 1(a), 
one spectrum would be interpolated using all the lines where N 2 is varied (blue circles), a second spectrum would be 
interpolated using all the 0 2 lines (black x’s), etc., repeated for all the species. Each of these spectra will be similar, 
but the n th spectrum will be more accurate in the vicinity of the resonances of the n th species since it is interpolated 
from the lines in the « th species direction (which have more closely spaced points). 

The first step in computing the 77 th spectrum described above is to interpolate along each line in the « th species 
coordinate, using quadratic interpolation, a spectrum at the desired concentration of that species. This results in a 
regular “reduced” grid of spectra at the node points of the remaining species, T nh where i=i'i,i 2 ,. . . ip.. .i^d where j=n 
is omitted; the reduced grid has A^-l dimensions and contains (X Ndl spectra, where O is the number of nodes in each 
coordinate direction. Since the reduced grid is regular, the 77 th spectrum is found using Lagrange polynomial 
interpolation as defined by the equation below: 
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The summation and product of a sequence in the above equation is omitted for j=n; y , is the composition being 
interpolated to and y iU are the compositions at the reduced grid points. Thus, for 4 resonant species and 3 nodes in 
each direction, as in the example of Fig. 1, the reduced grid is three-dimensional and there 27 terms 111 the 
summation; each summed term consists of the product of 27 terms that multiplies a theoretical spectrum value. The 
general formulation of this equation allows for efficient programming via “do” loops for arbitrarily large numbers of 
species. 

The final interpolated spectrum is the sum over the species of these spectra, weighted by the weight function 
distribution defined below: 


Z Nd 
n = 1 


WT 

n n 


The resultant spectrum is accurate for all resonances. 


C Residual Functions 

Experimental spectra are analyzed to determine composition and temperature by interpolating from the library of 
theoretical spectra and minimizing the residual functions defined below. S k is the square root of the experimental 

signal intensity (in square-root counts), 7 * is the square root of the normalized theoretical signal intensity, (Jjl is the 

variance of the experimental noise in the signal intensity (in counts), k is the detector pixel number (wavelength), N p 
is the number of pixels, and v is the number of degrees of freedom (N p minus the number of fit parameters). A is a 
scaling factor that rescales the theoretical spectrum so that it best fits the experimental spectrum, and is recomputed 
at each evaluation of the residual. The equation used to calculate X is slightly different for each residual; they are 
derived by differentiating the definition of the residual with respect to X, setting the result to zero, and solving 
analytically for X. 

Three residual functions are considered. The first residual is the mean square error of the signal intensity: 


Ri = 
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The resulting fit then is simply a least squares fit to the intensity. The second residual is the “reduced chi-square 
statistic”: 


V tt 

The resulting fit is then a weighted least squares fit to the intensity. Assuming that the noise model is correct, a 
value of » 1 indicates that the fit is poor and £ ~1 indicates a good fit. If £ <1 the model has too many free 
parameters and is fitting the noise. Snelling et al. 8 showed that both the random and bias errors in single-shot CARS 
measurements of temperature were dramatically reduced with weighted least-squares fitting compared to un- 
weighted. The third residual is the mean square error of the square root of the signal intensity and the resulting fit is 
simply a least squares fit to the square root of the mtensity. 

iy/ p 1 k=l 

The CARS signal intensity is roughly proportional to the square of the mole fraction of resonant species. Least 
square fitting of the square root of the signal intensity emphasizes the lower concentration resonant species and the 
non-resonant regions of the spectrum compared to least squares fitting the signal intensity. 

All previous single- and dual-pump CARS studies at NASA Langley (e.g.. Refs. 3, 4, 5) used least squares 
fitting to square-root intensity. (CARSFT has the options to least squares fit square-root intensity or a weighted 
square-root intensity.) The effect on the fitted parameters of the choice of residual function when fitting noisy CARS 
spectra will be studied in this paper. 

D Fitting Algorithm 

The fitting algorithm searches the library, interpolates spectra and minimizes the residual function. Fitted 
parameters are temperature, mole fraction of each resonant species, and a horizontal shift parameter that is applied to 
the data to correct small constant wavelength shifts in the data, A (cm -1 ). 

The algorithm includes a quick fitter to get an estimate of the solution. The quick fitter performs a direct search 
by computing the residual for every library spectrum and interpolatmg the composition and temperature that 
produces the minimum residual. This provides a good initial guess for the solution which is then refined by an 
iterative fitter. (This two-step process tends to avoid the possibility of the iterative fitter not converging to the 
global solution as it might if given an arbitrary initial guess.) 

The iterative fitter employs the Levenberg-Marquardt algorithm (LMA), as implemented by Lahey Fortran, 
which minimizes the function defined below: 

*=/<y)=2;{F t <y)} 2 

k = 1 

F is a function of the species mole fraction and pixel number, and can be identified by comparison with the 
definitions for the residuals. Rather than include temperature as a parameter of the fit, the algorithm is implemented 
to find the best fit at several fixed temperatures (includmg both library values and at temperatures intermediate 
between library values) and the minimum residual is computed at each temperature. The temperature that produces 
the overall minimum is interpolated using quadratic interpolation from the residuals at the fixed temperatures. The 
resulting algorithm was found to always converge to sensible parameter values that globally minimized the 
residuals. 

The horizontal shift is found in a separate iterative cycle applied both as part of the quick fitter and as part of the 
LMA fitter stages. The equation used to calculate A at each iteration is slightly different for each residual and is 
derived from the equation for the residual by setting T=T( v-A), where v is the wave number, differentiating this 
equation with respect to A, setting the result to zero, and then solving algebraically for A. The resulting equation 
contains dT/dA which is found by fitting a spline to T(v) assuming A= 0, which amounts to a local linear 
approximation to the slope. The equation is applied repeatedly until the accumulated shift ceases to change: the 
method typically converges withm 2-3 cycles. 
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III. Numerical Experiments 

Numerical studies were conducted to verify the accuracy of the new fitting algorithm. A set of spectra was 
generated numerically in which the concentrations of four resonant species (N 2 , 0 2 , H 2 , and C0 2 ) were varied at 
random and the spectra were fit to a sparse library using the algorithm. The Sandia code CARSFT (modified by 
Lucht 2 for dual-pump CARS and further modified by O’Byme 4 ) was used to calculate all the spectra. The spectral 
range of computation for this study was 1950 cm' 1 to 2400 cm' 1 and the second pump frequency offset was 650 cm' 1 . 
The first pump line width was assumed to be negligible, the second was 0.07 cm' 1 , and the instrument line width was 
2.0 cm' 1 . The Voigt line shape model was used for all species except H 2 where the Galatry model with narrowing 
parameter 0.047 was used. 9 

The library contained 9720 spectra: there were 9 temperatures from 800 K - 2400 K in 200 K increments, 
and 10 compositions for each species that were non-uniformly spaced with mole fractions in the range 0-0.21 
for 0 2 , 0-0.79 for N 2 , 0-0.5 for H 2 , and 0-0.25 for C0 2 . The library is sparse - a fully packed library with these 
dimensions would contain 90,000 spectra. 

Two hundred synthetic “experimental” spectra were generated using CARSFT that were uniformly and 
randomly distributed in the same temperature and mole fraction ranges as the library, but constrained so that 
the sum of the resonant mole fractions did not exceed 1 and N 2 mole fraction was greater than 0.3 for 
temperatures greater than 1600 K; this last condition is arbitrary but justified based on the fact that high 
temperature is not possible without combustion of fuel and air, leaving as one of the products N 2 . These 
spectra were rescaled by an arbitrary factor to be in the units of the detector analog-to-digital (A/D) converter 
or “counts” and noise was added to the spectra to simulate experimental noise. 



Intensity (counts) 


Figure 3. Model for standard deviation of signal noise. 

The noise was assumed to be Gaussian distributed with variance dependent on experimental signal intensity 
S 2 according to a model that included the dark noise of the detector, photon shot noise, and broadband dye 
laser mode noise; 10 the variance (in counts) is given by the equation below: 

(j\ = (D/N e ) 2 +S}/N e +{SlC bbdt f Equation 1 

D is the dark noise (12.6 electrons), N e is the number of electrons per count (1.05), and C^i is the intensity of 
the broadband dye laser noise as a fraction of intensity (0.115). These model coefficients were initially 
believed to simulate an existing CARS system at NASA Langley (see Section IV). Dark noise and broadband 
dye laser mode noise are expected to be Gaussian distributed. Shot noise obeys Poisson statistics that tend to 
Gaussian as the signal level increases (signal counts >10). Figure 3 shows the noise standard deviation plotted 
as a function of signal intensity, and also the individual contributions of dark noise, shot noise and dye laser 
noise. The ratio of total noise to signal does not go to zero but tends to a constant as the signal counts become 
very large; dye laser mode noise is dominant very much above 100 counts while the dark noise is dominant 
very much below this level. 
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The library took 75 hours to compute on a single processor of a PC, an average of 0.5 min/spectrum. 
Fitting the 200 experimental spectra to the library took about 150 minutes, 0.8 min/spectrum. If CARSFT in 
the fitting mode had been used, there would have been hundreds or thousands of spectrum evaluations to fit a 
single spectrum, taking many hours at best, if a converged fit could be found at all. 

Two typical experimental spectra and fits to them are shown in Figure 4: these are based on the same 
theoretical spectrum that has been rescaled to two different signal levels. The low signal-to-noise (S/N) case 
has 25 detector counts at the non-resonant background (NRB) baseline level and the high S/N case has 400 
counts at the NRB baseline level (note that the plot shows square root of signal level). Figure 4(a) shows the 
spectra for both cases and the fits to them using residual R 2 . Figure 4(b) shows the low S/N spectrum, the 
synthetic experimental spectrum prior to addition of noise, and the fits for each of the residuals. Note how fits 
using and R 2 are almost identical to each other and close to the no-noise spectrum, while the fit using R 2 is 
low where the signal is low and the noise highest. The low fit is expected since the noise distribution of the 
signal at a given pixel has mean equal to the no-noise intensity, but the noise distribution of the square root of 
the signal is skewed towards zero, and has a mean less than the square root of the no-noise mean. Obviously 
the best fit is one that recovers the no-noise spectrum since it will have the correct values of the fit parameters, 
and R { and R 2 gives superior results in this case. 
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Figure 4. Fits to synthetic experimental spectra with noise: (a) low and high signal-to-noise ratio (residual R 2 
minimized), (b) low signal-to-noise ratio with all three residuals fitted and prior to addition of noise. 


Table 2. Mean error and standard deviation of error in fitted temperature and composition fitting 200 noisy 
spectra, (a) No noise (baseline), (b) R x minimized, (c) R 2 minimized, (d) R 3 minimized. 




T 

mole fraction 

fit error 

res 

(K) 

n 2 

o 2 

h 2 

C0 2 


Ri 

1.3 

0.0001 

0.0001 

-0.0003 

0.0003 

mean 

r 3 

1.2 

0.0003 

0.0002 

- 0.0001 

0.0004 


Ri 

5.4 

0.0012 

0.0014 

0.0014 

0.0013 

st dev 

r 3 

5.1 

0.0008 

0.0010 

0.0007 

0.0009 


results 

NRB 

T 

mole fraction 

Ri 

(counts) 

(K) 

n 2 

o 2 

h 2 

C0 2 


400 

2 

0.0014 

-0.0007 

-0.0001 

-0.0015 

mean 

25 

5 

0.0027 

-0.0007 

-0.0003 

-0.0014 


400 

25 

0.006 

0.012 

0.006 

0.017 

st dev 

25 

35 

0.008 

0.017 

0.008 

0.020 


results 

NRB 

T 

mole fraction 


results 

NRB 

T 

mole fraction 

r 2 

(counts) 

(K) 

n 2 

o 2 

H 2 

C0 2 


r 3 

(counts) 

(K) 

n 2 

o 2 

h 2 

C0 2 


400 

-2 

0.0007 

- 0.0001 

0.0000 

- 0.0001 



400 

3 

0.0010 

0.0006 

-0.0005 

0.0005 

mean 

25 

3 

-0.0004 

-0.0003 

-0.0016 

-0.0004 


mean 

25 

46 

0.0239 

0.0126 

0.0046 

0.0157 


400 

13 

0.003 

0.003 

0.003 

0.003 



400 

16 

0.003 

0.003 

0.003 

0.003 

^ stdev 

25 

29 

0.006 

0.011 

0.005 

0.009 


st dev 

25 

54 

0.011 

0.021 

0.009 

0.021 
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Table 2 shows mean and standard deviation of the error (fit value minus true value) for the fitted 
parameters for 200 experimental spectra at different temperatures. Table 2(a) is a case with no noise (Ri are R 3 
are minimized). Tables 2(b,c,d) are cases with noise where R u R 2 are R$ are minimized respectively. The 
standard deviation of the errors in the no noise case, attributed to errors by the fitting algorithm, are ~5 K and 
less than -0.0014 in mole fractions, while the mean errors are negligible. These errors are small as compared 
to errors in fitting noisy spectra and thus the fitting algorithm introduces no significant errors. Comparing 
Tables 2(b,c,d) to each other it may be seen that residual R 2 produces consistently the best results, both in the 
mean and the standard deviation of the error in both temperature and composition. For example, for the low 
signal-to-noise case the mean and standard deviation of error in temperature are 3 K and 29 K respectively for 
residual R 2 , but 46 K and 54 K for R 2 . Results are similar for the mole fraction errors. Results for Ri though 
worse than R 2 are better than for R 2 . 


IV. Measurements in Hencken Burner 


Experimental broadband dual-pump CARS data acquired in the H 2 -air flat- flame produced by a Hencken burner 
were fitted using the algorithm. The Hencken burner flame has previously been carefully studied; 11 it was found to 
be very nearly adiabatic and the composition and temperature could be calculated using adiabatic, equilibrium 
assumptions. Burner flow rates were between 44 and 63 standard liters per minute while H 2 -air equivalence ratios 
were varied in the range -0. 1-3.5. 

CARS data were acquired using the portable system at NASA Langley 5 described most recently in the 
companion paper by Magnotti 12 (elliptical beam shaping was not used in the present experiment). For the present 
experiment laser beams pulsed at 20 Hz and -10 ns pulse length, one at 532 nm (injection seeded, doubled 
Nd:YAG), one at 552.8 nm (narrowband dye laser with line width -0.07 cm' 1 ), and one at 604.5 nm (broadband dye 
laser with full-width half-maximum line width of 9.6 nm), are focused and crossed in a planar BOXCARS phase 
matching configuration. The signal is separated by dichroic mirrors and focused at the entrance of a 1 meter 
spectrometer with a 2400 lines per inch grating. Spectra were normalized by a spectrum (averaged over -200 
pulses) acquired in argon (which is non-resonant) to remove the effect of the broadband dye laser line shape. 
Library spectra were generated using CARSFT; the Voigt lineshape model was used for N 2 and 0 2 and the Galatry 
model was used for H 2 with narrowing parameter of 0.047. 



(a) 



Intensity (counts) 


(b) 


Figure 5. CARS Hencken burner data at three equivalence ratios: (a) averages of 200 spectra, (b) noise levels 

compared to the noise model. 

The noise model of Equation 1 was employed in evaluating residual R 2 . The parameters of the noise model were 
estimated from experimental spectra. From about 200 experimental spectra in the argon, which is non-resonant, the 
ratio of the standard deviation of the noise to the signal level (excluding the wings of the spectrum where the dark 
noise is significant) was found to be -0.115; this is a good estimate of Cbbdi- Noise properties were also estimated 
from Hencken burner spectra at nominal equivalence ratios of 0.2, 1, and 2. The average of 200 spectra at each of 
these equivalence ratios is shown in Fig. 5(a); these spectra have not been normalized by the argon spectrum since 
noise depends on absolute not nonnalized signal level. For each equivalence ratio, the noise error in each pixel of 
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10 spectra was found by applying a small correction to make the total signal of each spectrum (summed over all the 
pixels) the same as for the average spectrum, to compensate for shot-to-shot random variations in the laser energy, 
and then subtracting the average spectra. The errors were binned according to the signal level of the pixel and the 
standard deviation of the error was found for each bin. The standard deviation of the error is plotted as a function of 
signal level in Fig. 5(b) at each equivalence ratio. Also shown is the noise model with D = 8 electrons and plotted 
both with C bbd f=0. 115 (labeled “model”) and CW=0.23 (labeled “Cbbdl x2”). The “model” line is thus a model of 
the noise in the non-resonant argon spectrum. The Hencken burner spectra data agree with the “model” line at lower 
signal levels, where the signal is generated by the non-resonant portion of the spectra, and then deviate to higher 
levels where the signal comes from resonant portions of the spectra. Noise levels seem to asymptote to the top 
model line where C bbdr 0.23. The reason the broad-band dye laser noise is higher in the resonant parts of the spectra 
is that none of these spectra are fully resolved by the spectrometer; thus, the interaction between the spectrally 
unresolved structures of the dye laser spectrum and the unresolved lines of the resonant spectrum produces greater 
noise than interaction with the featureless non-resonant spectrum. For the purpose of computing R 3 an intermediate 
value is used: C bbdr ~0.\6. Both theoretical and experimental spectra were multiplied by the reference argon 
spectrum to recover spectra in the original intensity units before computing R 2 : spectra in this form were used both 
directly and indirectly in computing the standard deviation of the noise error from Equation 1 . 

Problems were encountered in accurately fitting the isolated H 2 rotational lmes (the S(5), S(6), and S(9) shown in 
Fig. 5(a)). The apparent spectral location of these lines deviated a small amount from the location of the lines in the 
library spectra after fitting N 2 and 0 2 (error is less than the line width), even after correcting for small nonlinearities 
of the spectrometer with respect to wave number. Also the line shape was not always well represented by the library 
spectra. Since these are isolated rotational lines, at the resolution of our instrument the line shape contains no useful 
information. This problem was solved by applying the average level m a narrow range around the line center to 
each pixel in that range: this averaging was done to both the theory spectrum and the experimental spectrum prior to 
each residual evaluation. The result is that the algorithm fits only the area under the spectrum m the specified wave 
number range. An additional problem encountered with the H 2 rotational lines was that the temperature that 
produced the best fit to these lines was different than the calculated flame temperature. Several researchers report 
errors measuring temperature using H 2 CARS that are associated with line-shape modeling (e.g.. Ref. 11). If all 
three H 2 lines plus the N 2 and 0 2 resonant regions are included m the fit then the fitted temperature is some 
compromise between the temperature that best fits the H 2 lines and the temperature that best fits the N 2 and 0 2 
resonant regions; this compromise temperature depends strongly on which residual is minimized. Since the spectral 
modeling of N 2 and 0 2 has been widely validated for fitting for temperature (e.g., Ref. 1 1), and since S(5) is stronger 
than S(6) or S(9) over all our temperature range, only S(5) is included in the fitting; the other lines are omitted from 
the residual by excluding small wave number regions each side of each line center. The result is that H? has no effect 
on the fitted temperature but the algorithm fits for H 2 mole fraction. Another problem is that none of these H 2 
rotational lmes are significantly populated when room temperature is approached (between 300K and 500 K). At 
these temperatures the method becomes ill-posed for composition and small noise errors in the spectrum are fitted as 
large mole fractions of H 2 . Tedder et al. 7 have proposed using the S(3) line which is strong at both room and flame 
temperatures. 

Figure 6 (a-d) shows the mean and standard deviation of the fitted temperature and the species mole fractions 
from 200 spectra at several equivalence ratios plotted as a function of equivalence ratio. Results are shown for fits 
minimizing each type of residual. For the purposes of the plots, the standard deviations are multiplied by 10. Also 
shown is the calculated temperature and mole fraction distribution based on flow rates to the burner; uncertainty 
bands on the calculations are associated with uncertainties in the flow rates. The error between theoretical and mean 
measured temperature is typically less than about ±50 K and the error in mean mole fractions is less than about 
±0.02. Mole fraction of 0 2 is measured high, although the results of minimizing residual R 2 appear to be the best. 
Errors in N 2 and Fl 2 are similar for all three residuals. 

The standard deviation in fitted temperatures and mole fractions are consistently lowest for residual R 2 , followed 
by i? 3 and then R\. In some cases the trend is very pronounced; in the case of temperature or 0 2 mole fraction the 
reduction from Ri to R 2 is a factor of one half or less. This is an important result for the use of CARS to measure 
turbulent fluctuations in flames since instrument random error limits the minimum turbulent fluctuation level that 
can be resolved. These results are consistent with the results of the numerical experiments in Section III in that 
minimizing the chi-square residual, R 2 , reduced random error (noise) and bias errors in some fitted parameters. The 
standard deviations for this residual are -50K for temperature, -0.005-0.0 1 for 0 2 , -0.01-0.015 for N 2 , and -0.01- 
0.03 for H 2 . The Hencken burner flame may have some real fluctuations in temperature and composition, so the 
random errors in this CARS method may be less than suggested by these results. 
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Figure 6. Fits to CARS Hencken burner data; mean and standard deviation of fitted parameters as a function 
of equivalence ratios. Standard deviation is multiplied by 10. (a) Temperature, (b) mole fraction 0 2 , (c) mole 

fraction N 2 , (d) mole fraction H 2 . 


V. Conclusions 

This paper has described the development of an algorithm for fitting experimental CARS spectra containing 
multiple resonant species to theory. The algorithm employs a library of previously calculated theoretical spectra and 
thus is much faster than methods that compute new spectra at each iteration of the solver. Library spectra are 
assumed calculated by some other algorithm. The library structure is orthogonal in a non-resonant susceptibility 
weighted mole fraction coordinate and is sparse; its size thus grows much more slowly with number of species than 
a conventional fully-packed library. The method of interpolation from the library is discussed; several residual 
functions minimized by the fitter were considered, and the fitting method was described. Synthetic “experimental” 
dual-pump CARS spectra containing 4 resonant species (N 2 , 0 2 , H 2 and C0 2 ), both with noise and without it, were 
fitted by the algorithm and the method was demonstrated to be both robust and accurate (there were negligible errors 
in fitting noise-free spectra). Experimental dual-pump CARS from a Hencken flat-flame burner reacting H 2 and air 
were also fitted. In both studies, weighted least square fitting of signal intensity, as opposed to least squares fitting 
signal intensity or square-root signal intensity, was shown to produce the least random error and the same or less 
bias error in the fitted parameters. Standard deviation of the error (random error) was reduced by factor of one half 
or more in some parameters by weighted least squares. Use of weighted least squares is thus particularly important 
in turbulence studies where minimum resolved turbulence levels are limited by random error. Some experimental 
issues associated with fitting isolated rotational lines of H 2 were also resolved. 
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